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Abstract. Boundary problem for linear partial differential algebraic equations system with multiple 
characteristic curves is considered. It is supposed that matrix-functions pencil of the system under 
consideration is smoothly equivalent to special canonical form. For this problem, with the help of the 
spline collocation method, a difference scheme of arbitrary degree of approximation with respect to 
each variable is constructed. Sufficient conditions for its absolute stability is found. 



1. Introducton 

When modeling some processes of hydrodynamics, gas dynamics, atmospheric physics, plasma 
physics etc., it appears systems with identically degenerate matrix- functions in its domain of defi- 
nition at all higher partial derivatives [I]-[E]- Such equations are known in the literature as partial 
differential algebraic equations, equations not resolved with respect to higher derivatives, degener- 
ate systems of partial differential equations or Sobolev's equations. There exists a number of lines 
of research of these equations. One of these lines is based on investigation of canonical structures 
of matrix pencils [3J, [9], [HI]- [13]. At the moment in the literature, it is well studied the issues of 
the existence and numerical solution of partial differential algebraic equations with constant matrix 
coefficients [3J, [I], [8], [H]. It is explained by the fact that canonical structures of constant matrix 
pencils are well studied. Insufficiently studied global properties of matrix-functions pencils impede the 
investigation of the systems with variable matrix coefficients. In the works [10]- [13] we considered the 
systems with simple characteristic curves whose matrix-functions pencils satisfy the criterion known 
as "rank-degree" or double "rank-degree" one. Matrices in these pencils depend only on two variables 
|15| . For such systems we constructed two-layer and three-layer implicit difference schemes of a first 
and second orders of approximation |10|-|13|. This accuracy is not always enough, for example, in the 
case when the Lipschitz constant and a domain of definition are quite large. In this case, it requires 
difference schemes of higher orders of approximation. In the work [16j . we studied matrix- functions 
pencils depending on many variables. In a result we obtained sufficient conditions of a smoothly equiv- 
alence of these pencils to canonical structure similar to Kronecker form of a regular pencil of constant 
matrices. This structure is in fact a generalization of a canonical form of the pencil satisfying to the 
"rank-degree" criterion. This structure allows to investigate linear degenerate hyperbolic systems with 
arbitrary number of independent variables with multiple characteristic curves. 

In this paper we consider linear partial differential algebraic equations system of hyperbolic type 
with the pencil smoothly equivalent to its canonical form similar to the Kronecker form [16]. With 
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the help of spline collocation method, the foundations of which are set out in [T7J, we construct a 
high-performance implicit difference scheme of higher order of approximation and then we prove its 
stability. 

The rest of the paper is organized as follows. In the next section, we give a statement of the 
problem. In section [3l we perform an approximation of unknown function on the uniform grid by 
a spline of degree mi and m.2 with respect to each variable, respectively. Then we write down our 
difference scheme. This approach for constructing of difference scheme was used in [13]. In the work 
|13j . for an approximation of unknown function we used the bi-cubic spline with defect. Accuracy of 
the calculations in this case was low. It is worth remarking that in the present work the corresponding 
results significantly improved. In section^ we give some notations and auxiliary propositions necessary 
to justify correctness of difference scheme are given. The section [5] is designed to cast the difference 
scheme to some canonical form. In this form one can easily analyze the spectrum of matrix coefficients. 
It is shown that for any values of steps, the spectrum is entirely contained in the circle of unit radius. 
In section [6] we prove an absolute stability property of our difference scheme. Finally, in section [71 for 
two test examples we show results of numerical experiments. 

2. Statement of the problem 

Let us consider a boundary problem for the linear partial differential equations system 

A(x, t)d t u + B(x, t)d x u + C(x, t)u = f, (2.1) 
u(x ,t) = ip(t), u(x,t ) = (t>(x), (2.2) 

where A(x,t), B(x,t) and C(x,t) are n x n matrices whose elements are supposed to depend on 
iGl 1 and t € R . In the following we suppose that (x, t) G U = [xq; X] x [to; T] and U is therefore a 
domain of definition of all functions under consideration. It is supposed that the elements of matrices 
A(x,t), B(x,t), C(x,t) and free term f(x, t) belong to C 2 (U). The vectors tp(t) and 4>(x) are supposed 
to be some given vector- functions of its arguments whose elements belong to C 2 ([to, t]) and C 2 ([xo, X]), 
respectively. 

Let us suppose now that 

det A(x,t) = and det B(x,t) = V (x,t) € U. (2.3) 

The system (|2.ip with the condition (|2.3p is therefore partial differential algebraic equations one. Its 
investigation is closely related to analysis of global properties of the pencil A(x,t) + XB(x,t). In 
the work |16j . we obtained the conditions of smooth equivalence of matrix- functions pencil to its 
canonical form similar to canonical structure of regular pencil of constant matrix, in this connection, 
let us remember the definition of smooth equivalency of matrix-functions pencils. 

Definition 2.1. |16| Two n x n matrices pencils A(x, t) + XB(x, t) and A(x, t) + XB(x, t) with elements 
belonging to C S (U), where A is a some parameter are called s-smoothly equivalent if there exist square 
matrices P(x,t) and Q(x,t) which do not depend on A and satisfy following conditions: 

(1) the elements of matrix P(x,t) and Q(x,t) belong to C S (U); 

(2) V(x,t) £ U there exist P _1 (x,t) and <5 _1 (x,t); 
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(3) the relation P(x,t)(A(x,t) + XB(x,t))Q(x,t) = A(x,t) + XB(x,t) V(x,t) G U holds. 

It is useful following theorem |16| . 
Theorem 2.2. [16] Let the following conditions be fulfilled: 

(1) all roots of the characteristic polynomial det(A(x,t) + \B(x,t)) are real and have a constant 
multiplicity in a domain of definition U; 

(2) ranks of matrices A(x,t) and B(x,t) are constant at each point of a domain U and less than 
n. 

Then the pencil A(x,t) + X B(x,t) is smoothly equivalent to the canonical one 

diag{£ d , M(x, t),E p } + A diag{ J(x, t),E h N(x, t)}, (2.4) 

where Ed is an identity matrix of an order d; M(x,t) and N(x,t) are an upper (right) triangu- 
lar blocks with zero diagonal of orders I and p, respectively; Oi is a zero square matrix of order I; 
J(x,t) = diag{ J\(x, t), J2(x,t), . . . , Jk(x,t)}, where Ji(x,t), for i = 1, . . . , k are nonsingular matrices 

k 

of orders d%, respectively; d= ^ d u ; each block Ji has a unique eigenvalue —1/Xi(x,t) in the domain 

of definition U; Xi(x,t) are eigenvalues of characteristic polynomial det(A(x,t) + XB(x,t)), different 
from zero in the domain of definition U ; p = n — d — I. 

Remark 2.3. Let conditions of the theorem 12.21 are valid. Then if one requires 

rank(5(x, t)) = deg (det(A(x, t) + XB(x, t))) 

or 

rank(A(x, t)) = deg (a det(A(x, t) + B(x, i))J , 

then the pencil A(x, t) + XB(x, t) is s-smoothly equivalent to the pencil (|2.4p . in which N(x, t) = Q[ or 
M{x, t) = Op, respectively. In this case one says that the pencil A(x, t) + XB(x, t) satisfy the criterion 
"rank-degree". A structure of such a pencil was investigated in [15] . 

Let us suppose that matrix-functions pencil A(x,t) + XB(x,t) of the system (|2.ip satisfy the con- 
ditions of theorem 12.21 and remark 12.31 In the next section, we proceed to construct our difference 
scheme. 

3. DIFFERENCE SCHEME 

To construct the difference scheme , we perform the partition of the domain U by the lines x = x% , 
where xi = xq + ih and t = tj, where tj = to + jr for i = 0, . . . , n\ and j = 0, . . . , n-z. In a result, 
we obtain a uniform grid XJ/\ with the steps h and r with respect to the space and time variable, 
respectively. Clearly, h = (X — xq)/ui and r = (T — io)/ n 2> where < r < To and < h < ho. The 
points with coordinates (x{,tj) are referred to as the nodes of the grid while the lines x = x, t and 
t = tj are called its layers. In each domain Uij = [xi,Xi + m\h\ x [tj,tj + m2r] C Ua , where m\ < n\ 
and m2 < n,2, we seek the approximation of the solution u(x,t) of the problem (|2.ip . (|2.2p in the 
form of a polynomial L™ 1,m2 (x,t) with indeterminate coefficients of degree mi and m2 with respect 
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to variables x and t, respectively. We require that the values of this polynomial L™ 1,m2 (x, t) in the 
nodes (xi + l±h, tj + 1 2 t) for l\ = 0, . . . , mi and l 2 = 0, . . . , m 2 of a grid domain U{ coincide with the 
values of the desired solution u(x, t) at these nodes. In order our approximation will be continuous, 
we require that the relations 

LT}'™ 2 (xi + hh,t 3 ) = + hh,tj), L™(x o ,t 3 = ^%), 



'•j 



T m 1 ,m,2 
Hi 



Xi,tj + l 2 r) 



"^j^ 2 (xi,tj + 1 2 t), L^' m2 (xi,t ) = 4>(xi). 
hold at horizontal and vertical line x = x% and t = tj, respectively. Applying differenceless formulas 
of numerical differentiation for equidistant nodes ([IS], p. 161), to derivatives dtu(x,t) and d x ii(x,t) 
at layers x = x% and t = tj, we obtain 



1 



d x u(xi + hh,tj) = d x Lfj' m2 { Xi + hh,tj) + ei(/i mi ), h 
d t u{ Xi ,tj + l 2 T) = d t L™y m \xi,tj + 1 2 t) + e 2 (r m2 ), Z 2 = 1 



,...,mi 
...,m 2 , 



(3.1) 
(3.2) 



where 



, mi 

d x L™l' m2 (xi + ZiMj) = - ^2 %,h u ( x i + hh,tj), 

mi 

d t L™J ,m2 (xi,tj + 1 2 t) = - ^ ^ hth u(xi,tj + Z 3 t), 

, Xj < Ci < a;i+i, 



ei(/* mi ) 



«3=0 



(mi + 1)! 



v=0 
rn 2 



T m2am2+lf T . /• \ /fljm 2 +l j 

e 2 (r mi ) = T ^ 7 ^'C2)/Of _d_ yj (s _ , 

' 1 as ■ L1 



(m 2 + 1)! 



u=0 



S = ll 



S = lo 



, tj < ( 2 < tj + l. 



The coefficients "fi lt i 3 and 7i 2 ,i 3 are calculated with the help of following formulas ([19J, p. 161): 



Huh = H(m,s,l 3 )\ 



m=m\ , 



C l * d 



(m, s, i 3 ) = {-l) m+h \ J[{s- u)/{s - I 



(3.3) 



\u=0 



where stand for the binomial coefficients. Writing the system (|2.ip at the nodes (xj + l\h, tj + Z 2 r) 
for li = 1, . . . , mi and l 2 = 1, . . . , m 2 of the domain Uij and substituting into it the values of desired 
function u(x, + l\h,tj + 1 2 t) and approximation of its derivatives (|3.ip and (|3.2p in these points, we 
obtain following difference scheme: 

m 2 . mi 

h=o h=o 
+Ci+i lt j+i 2 Ui + i 1 j + i 2 = f i+ i lt j + i 2 , 

uo,j=ipj, u ifi = (t>i, i = 0, . . . ,n 1 - 1, i = 0, . . . ,n 2 - 1, (3.4) 



where 



A i+h j + i 2 = A(xi + hh, tj + Z 2 r), B i+ i uj+ i 2 = B{xi + /i/i, tj + / 2 r), 
Cj +Zlij+Z2 = C(xi + + Z 2 r), f i+h>j+ i 2 = f(xi + + 1 2 t), 

u i+hJ+h = u ( x i + hh,tj + / 2 r), <fo = 4>(xi),if)j = tp(tj). 
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The difference scheme (|3.4[) at each node of the grid is given by a system of linear algebraic 
equations of order n = mim 2 n with unknown vector 

^i+lj + l — (^i+lj'+l) • • • i ^i+l,J+m2> ^i+2,j'+lj • • • j ^i+2,j+m 2 : • • • j "t+mu'+li • • • j ^£+mi +m2 ) • 

In what follows to avoid a confusion in notations, we denote the solution of the system (|3.4p in the 
nodes (xj+i, ij+i) of the grid ?7a by while the values of unknown function u(x, t) in the same 

nodes we denote as u(xi + i, = lij+ij+i. Then the system (|3.4p has the following form: 

m 2 ^ mi 

= fi+h,j+l 2 ~ -Ai+l u j+l 2 Jl 2 ,0 v i+h,j ~ -^Bi+h,j+l2%fi v i,j+l2i ( 3 - 5 ) 
vo,j = ipj, v ifl = 4>i, i = 0,...,m-l, i = 0, . . . ,n 2 - 1, li = l,...,mi, l 2 = l,...,m 2 , 
where 

= ('Wj+lj'+l) • • • > v i+l,j+m 2 , v i+2,j+l, ■ ■ ■ ; j+m 2 > • • • j^i+mu'+lj • • • ; ^j+mi ,j+m 2 ) T • 

Let us remark that the difference scheme (|3,5p in fact represents the whole set of implicit difference 
schemes. Setting different values of the orders mi and m 2 of approximating polynomials L™ 1 '" 12 (x, t), 
we get according to (|3.ip and (|3.2p difference schemes with corresponding approximation orders 
0{h mi ) + 0(r m2 ). In addition, the scheme (|3.5p can be both two-layer and multi-layer depending 
on given information. 

In our case we have in mind two-layer difference scheme, that is, we suppose that the values of 
unknown function are given only at one left layer x = Xi and on one lower layer t = tj of Ua (at 
the nodes of these layers the values of grid function were defined in the previous step), but in actual 
calculations are used mi + 1 vertical and m 2 + 1 horizontal layers, as is shown on the template (Fig. 
1). The movement on the grid is carried along its layers and within each layer by steps. Therefore we 
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FIGURE 1. The grid Ua and (mi + l)(m 2 + l)-point template. 

have constructed the difference scheme (|3.5p or more exactly a whole set of difference schemes with 
different orders of approximation. To proceed, we must to find out the conditions of solvability for 
the difference scheme (|3.5p for any data and for sufficiently small grid steps and finally to prove its 
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stability following the definition from ([18]. p. 12). It is worth remarking that the difference scheme 
(|3.5p was written in [20] but without proof of its stability. 

To estimate a solution of the problem (|3.5p . in the space C(Ua) consisting of n-dimensional vector- 
functions v(xi, tj) = v 2 (xi, tj), . . . , v n (xi, tj)) T let us use the following norm: 

\\v(xi, tj)\\cai A ) = max 11^(^)^)11) where \\v(xi,tj)\\ = max \v u (xi, tj)\. 

(xi,tj)£UA u=l,...,n 

It is obvious that it is consisted with the norm of n-dimensional vector- function v(x,t) £ C(U) [23] 

IKM)llc(tO = max{||u(x,t)|| V(x,i) G U}. 

In the following, we need in some auxiliary propositions which we formulate and prove in the next 
section. 

4. Preliminaries 

Let us spend some lines to fix auxiliary notations used throughout the paper. Denote r/h = r and 
define diagonal matrices 

7mi = dia g{7i,o,72,o, • • • ,7mi,o}, 7 m2 = diag{7i j0 ,72,o, • • • ,7m 2 ,o}, 

lmi = (%j)T\ 7m 2 = (li,j)T 2 (4.1) 

whose elements are given by (|3.3p . Remark, that j mi and j m2 coincides if m\ = m-2- Also we define 
following block diagonal matrices: 

mathcalA = diag{A 1 ,A 2 ,A 3 }, Si+ij+i = diag{S^ +1 S 2 , !B 3 }, (4.2) 



where 



= E mi ® 7m 2 ® + r7 mi (8) £ m2 ® Ji+ij+i, 
^ 2 = ^7 mi <8 £' m2i , ft 3 = £J mi (g> 7 m2 (8) 
.A 1 = E mi 7 ^ 2 (8 E d , A 2 = O mim2h A 3 = E ni ® 7 ^ 2 ® £ p , 
= 7m! ® ^m 2 <E> Ji+ij+i, = 7mx ® E m2 ® E h !B 3 = O mi 



m 2 p> 



Lemma 4.1. Lei ^ j , and J+1 /or si = 1, . . . , mi, s 2 = 1, . . . , m 2 and s 3 = 1, . . . , k are 

eigenvalues of the matrices j mi , 7m 2 and Jj+ij+i, respectively. If at each node of the grid U/\ the 
inequalities 

r^ 1 ^-C 2 , ^8i,8 2 and 8 3 , (4.3) 

be fulfilled, then the matrix f)j + ij + i in fl^.£[ ) ^ s nonsingular on the grid Ua- 

Proof. To prove the nonsingularity of the matrix on Ua under condition of this lemma, 

it is sufficient to prove this property for its blocks. Since the matrices j mi and 7 m2 are nonsingular 
then the blocks and J7 3 also share this property. Let us consider now the block Let "?% 

and CP2 are constant matrices and ^(a;, t) is a matrix- function, such that 

PrTmiOT 1 = V. ^Tm^j 1 = 7m 2 ) mX,t)J(x,tp^(x,t) = J(x,t), 
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where 7 mi , 7 m2 and J(x,t) are normal Jordan forms of the matrices 7 mi , 7 m2 and J(x,t), respectively. 
Remark that an existence of the matrix ^(x, t) is provided by the theorem in [22]. Make up 7(x, t) = 
"?i <8> J > 2 <8> y^x, i) and multiply the matrix - +1 on the right and on the left by matrices ^Pi+i.j+i 
and respectively. In a result, we obtain a normal Jordan form of the block It is 

obvious, that provided (|4.3[) there are no zeros at the main diagonals of Jordan cells of normal form 
of the matrix Therefore the lemma is proved. 

The following lemma is particularly important because it gives characteristic properties of splines 
and allows to transform the difference scheme (|3.5|) to its canonical form in which one easily sees 
spectra of matrix coefficients. 

Lemma 4.2. Let xq £ M. d be an arbitrary vector and yo = (xq, xq, . . . , £o) T £ M. md . Let a be an 
arbitrary parameter and 3 be an arbitrary constant d x d matrix with corresponding eigenvalues Q for 
I = 1, . . . , d satisfying Q / — /a, V /, s ; where /or s = 1, . . . , m are eigenvalues of m x m 
matrix 7 m = ('Ji,j)i l with elements defined through \3. 3\) . Then 

(E md + a(7" 1 ® a))" 1 y = -diag {exp(-aa), exp(-2a0), . . . , exp(-maa)} y + 0(r m ). (4.4) 



Proof. Let us consider Cauchy problem for homogeneous system of ordinary differential equations 

ct 

x(t) = 3x(t), x(t ) = x , t£/=[0,mr]. (4.5) 

r 

The solution of (14.51). as is known, is 



x(t) = exp (--3*) (4-6) 

Divide the segment / into (m — 1) equal parts by the points tj = jr for j = 0, ... ,m. Making use 
of differenceless formulas of numerical differentiation (|3.2p . we write down the approximation for the 
derivative of the vector- function x(t) at the points tj for j = 1, . . . , m 

m 

Tx(jT) = Y J W( t i) + °(T m ). (4.7) 

1=0 

Rewrite the relations (j4.7j) in the matrix form 

Y = ( 7 ° m E d ) y + ( 7 m ® E d ) y + 0(r m+1 ), (4.8) 

where Y and y are vectors of the form Y = (x(h),x(2h), . . . , x(mh)) T and y = (x(h),x(2h), . . . , x(mh)) T , 
respectively. With the help of (|4.5p we find the values of the derivative of unknown function at the 
points tj of the segment /. Substituting these values into (|4.8|) . we get 

- a(£ m ® a)y = (7m ® ^) 2/0 + (7m ® S d ) y + 0(r m+1 ). (4.9) 

Resolving (|4.9p in favor of y gives 

y = - (a(E m ®3) +lm ® Ed)- 1 ( 7m 8) E d ) y + 0(r m+1 ). (4.10) 

In the work |24j . we proved the identity 

(l™ 1 l Q m®E d )e md = -e md + 0{T m ). (4.11) 
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Here e m( j stands for a md-dimensional vector, each element of which ia a unit. Remark, that relation 
(|4.1ip is in fact a particular case of the following one: 

{imlm ® Ed) Vo = -V0 + 0(r m ). (4.12) 

Remark, that it is proved in the same way as ()4.1ip . but instead of the interpolated vector- function 
g{t) = (1 + t m , 1 + t m , . . . , 1 + t m ) T , we take here g{t) = x - e d + g(t). Making use of (kTT2]) . we 
transform (j4.10p to get 

V = ~ {E md + a( 7 " 1 ® a))" 1 yo + 0{r m ). (4.13) 

Finally, substituting into (|4.13p the values of unknown function (|4.6p at the nodal points of /, we 
obtain the relation (|4.4p . Therefore the lemma is proved. 

5. Transformation of difference scheme under consideration 



It is quite difficult to investigate the difference scheme (|3.5p without preliminary transformation of 
it to convenient form. Thus, let us first to cast it to special canonical form. As a result we write down 
this canonical form in the end of this section. 

We suppose that the pencil A(x,t) + XB(x,t) of the system (|2.ip satisfies the conditions of the 
theorem 12.21 Remember that this means that one can find a pair of matrix-functions P(x, t) and 
Q(x,t) with properties of definition 12. 1\ which serve for transformation of our pencil to the canonical 
form (|2.4p . Prepare following square matrix of order n: 

P = diag {-Pj+ij+i, . . . , Pi+ij+ m2 , Pi+2,j+li ■ ■ ■ i Pi+2,j+m 2 i ■ ■ ■ > Pi+mi,j+li ■ ■ ■ i Pi+mi ,j +m 2 } • 

Multiply the left- and right-hand sides of the system f)3.5[) on the left by the matrix rP and per- 
form a change of variable Vij = QijWij, where Wij is some unknown n-dimensional vector-function, 
calculated in the node (i, j). In a result we obtain the following difference scheme: 

( mi mi 

Pi+h,j+hi Ai+hd+h lh,hQi+h,j+h w i+h,j+h + r Pi+h,j+h %,l3Qi+h,j+h w i+h,j+h 

+TC i+ i 1 j + i 2 Q i+ i 1 j + i 2 w i+ i 1 j + i : }j = P i+ i lt j + i 2 ^Tfi + i u j + i 2 - rAi +h j + i 2 ^i 2fi Qi +h jW %+h j 

- r Bi+ii,j+i2%,oQi,j+i2 w i,j+h^ji 

for i = 0, . . . ,m — 1, i = 0, . . . , ri2 — 1, h = fn\ and I2 = 1, • • • , vcii with unknown vector 

"^j+lj'+l = (^i+lj'+l) • • • i w i+l,j+m 2 i w i+2,j+l, ■ ■ ■ ,Wi+2,j+m 2 i ■ ■ ■ j w i+mi,j+l, ■ ■ ■ , ^i+mi ,j+m 2 ) T • 

With Taylor's formula we can represent the matrix Qi + i ± j + i 2 in the following form: 

Qi+luj+h = Q i+ i 1)j+ i 2 + <?l + i uj+ i 2 Xih + crf + i uj+ i 2 X2T, 
for li = 1, . . . , mi an d I2 = 1 • • • ; m2, where 

a l+h,j+h - 9xQ{xi+h + °K tj+h) and a i+h,j+i2 = d T Q(x i+h ,t j+ i 2 + 0t). 
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Here 6 is some parameter, by assumption, obeying the condition < 6 < 1 and 1^ are suitable numbers 
for which one has Xk = \h — h\i while < Xk < m^. Let us denote Cij = PijCijQij and rewrite 
system (|5,ip . taking into account canonical form (|2.4p of the pencil A(x,t) + XB(x,t). In a result, we 
obtain 

^i+lj+l^i+lj+l = <7i+lJ+l> (5-2) 
^i+lj+l = + /iCTl + T(72, 

ft+lj+i = - (A + Ta 3 )w i+ ij - r(B i+ ij + i + To- 4 )-u)jj + i, 

where Oj+i ,i+i, .A and ^B^+ij+i are n x h matrices of the form 

Oi+ij+i = # mi <g 7 m2 <g diag-t-Ed, ; , £p} + rj mi (g) £ m2 <g diag{ J i+ i J+ i, £7, P }, 
yi = E mi <g 7^ 2 <g diag{£ rf , 0/, E p }, l B i+1 j +1 = j mi (g E m2 <g diag{ Jj+i^+i, E h P }; 

and 

+ l = TPfi+l,j + l, 

fi+l,j+l = (/i+lj+lj • • • > fi+l,j+m2i fi+2,j+li ■ ■ ■ i fi+2,j+m,2i ■ ■ - > /i+mij+lj • • • > /i+mij+rfi2 ) i 

Wj+ij = (e m2 ® Wi+i,i> e m 2 ® w i+2,i, • • • , e m2 ® Wj +mi j) T , 

= e mi <g (Wij+l, • • • , Wi,j+m 2 ) T ■ 

Here at are bounded in the domain of definition matrices obtained as a result of expanding Qi + i lt j + i 2 
and Jj + / 1J+ / 2 by Taylor's formula. 

Let us split each block component Wi + ij + ± of the vector u>i + ij + i into three ones 

= ( w i+l,j+li w i+l,j+li w i+l,j+l) T 

of the size d, I and p, respectively. Prepare then following n x h permutation matrix 

T = colon(Ti, T 2 ,T 3 ), 

T 1 = E mi ®E m z®{E A Q Axl Q dxp ), T 2 = E mi ®E m2 ®(Q lxd E l Q l><p ), (5.3) 
T3 = E mi (g E m2 (g (Opxd OpxZ -Ep)) 
where the blocks T\, T2 and T3 are supposed to have the sizes mxm^d x n, mxm^l x n and m\m2P x n, 
respectively. It is easy to prove that T is orthogonal matrix. The matrix T perform the permutation 
of the components of the vector Wi + \j + \ so that 

_ rp - _ / -1 -2 -3 \T 

Zi+lj+l = J ^i+ij+i — ^i+lj+l, , 

where 

^i+lj'+l = ( w i+l,j+li • • • ' w i+l,j+ni2i W i+2,j+li • • • ) w i+2,j+m 2 ' '"J w i+m 1 ,j+l^ ■ ■ ■ > w i+m 1 ,j+m 2 ) ' 

One sees that the vectors +1 and have the sizes mim2d, m\vri2l and mirri2P, 

respectively. With the help of T we also reshuffle the elements of f2j+ij+i, A and Bj+ij+i in the 
following way: 

rOj+i^+iT = Oj+ij+i, T.AT T = A, T"Bi + ij + iT T = Bi+ij+i, (5.4) 
where Jlj+ij+i, ^4 and were defined in (|4.2p . 
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Multiply the system (|5.2p on the left by the matrix T and write the unknown vector Wi+ij+i as 
Wi + ij + i = T Zi+ij+x- Then, taking into account f)5.3j) and (|5.4|) . we obtain from f)5.2j) the following 
system: 

^i+ij+i^i+ij+i = rgi+ij+i - (A + ^3)^+1^- - r(S i+ ij + i + ra^Zij+i, (5.5) 

where 

^i+lj+l = ^i+lj+l +ha\ + T(T2,gi+l,j+l = Tfi + ij + i 

and <7fc = T<JkT T . 

Next we suppose that in each node of the grid [/a inequalities (|4.3p from the lemma 14.11 are valid 
and r being a ratio of two steps is a constant. In virtue of the lemma 14.11 the matrix Oj+ij+i is 
nonsingular on the grid Thus, the matrix f^+i .j+i in (|5.2p can be presented in the following form: 

where Mj + i J+ i = J5jj + J+1 (/i<7i + TCT2). Since, by assumption, the elements of matrix coefficients 
of the system (|2.ip belong to C 2 (U), then in virtue of the theorem 1 in the work [22], the elements 
of J(x,t) also belong to C 2 (U). This means that there is such a constant 3C, for which the condition 
imij+ill < 3C is fulfilled. 

In virtue of the theorem from ([21]. p. 195) there exists the inverse matrix It can be 

presented in the form - +1 = E n + Mj+ij+i/i, where 

00 

M i+W = J^(-l) 

!/=l 

Therefore we see that the matrix fij^ j +1 exists. It is worth remarking, that therefore we proved the 
solvability of the difference scheme (|3.5p in each node of the grid {/a for sufficiently small steps h and 
r. 

To proceed a transformation of the difference scheme (|3.5p . we multiply ([5.5]) on the left by the 
matrix j+1 to get 

= - + ei+i,j+i(h,T))z i+ ij - r(K i+ ij + i + e i+ ij + i(/i, r))zij + i, (5.6) 

where is the n-dimensional vector given by <?»+ij+i = riW^. ljJ . +1 f27 + 1 j+iflH-lJ+i I 

-^i+ij+i = -Ki+ij+i = ^i+ij+i^^+ij+i; 

and 

h+i,j+i(h,r) = hM i+1J+1 nr^ lj+1 A + TQr^ lj+1 d 3 + ZirMj+ij+iQ^ j +1 cr 3 , 
e i+ ij +1 (/i, r ) = ^M i+1 j + i0^. 1)i+1 S i+1 j + i + ril^ 1 ljj+1 a- 4 + /it M i+ i J - + iil l ^ 1 1J+1 a- 4 
being h x n matrices. 

Taking into account (|4.2p . the system (|5.6p is decomposed into three subsystems 



k Tpk „-,k 
t+lj+l — r i+lJ+l w i-\ 


-1J 


~~ r -^'i+i,i+i u; « fe ,i+i 




3 




3 










(5.7) 






!=1 





fr 1 



(era + r<7 2 ) 
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Here 



TTlfc 



o k 

s 'j+i,i+i 



1 a k K k 



Si i+l,j+l 



1 u 
rr>K 



while e^j-i J+1 (/i, t) anc ^ ^i+lj+i r ) are matrix blocks of 

e mj - + i(/i,r) = (ef+ij+iC^T)) and e m , i+1 (7t, r) = (e, fc + lji+1 (/i, r)) i , 

respectively. The sizes of these blocks correspond to the decomposition of Zi+\j+\. Finally, free terms 
9i+i,j+l are blocks of the vector 

9i+i,j+i = G?i+ij+i> £ff+ij+i> 5 4 3 +ij+i) T 

of the sizes mim^d, m\mil and m\m2P, respectively, and the vectors w k +1 j and w k j +1 have the 
following form 

Using (|4.12p . we transform the product F^ +1 j +l wj +1 j in the first equation of the system (|5.7p to 
get 

^VlJ+l^i+lJ = (-^1 ® 7m 2 ® + ^7mi ® ® Ji+lj+l)" 1 (^1 ® 7m 2 ® #d) 

= {E mim2 d + r7 mi <8 7^ (8) Ji+ij+i) 1 + 0(r m2 ). 

Let -R be a mi x mi matrix transforming j mi to the normal Jordan form. Remark that the matrix 
j mi is simple for any value of mi, that is, 

Rj mi RT 1 = f mi , where f mi = diag , , • • • , } • 

Then 

i^+lj+l^i+ij = ® E m 2 d) R(R® E ni2d ) w} +1>j + 0(r m2 ), (5.8) 
with R = diag{-Rii,-R 2 2, . . . ,R mimi }, where 

R ss = (E m2d + r^ mi (-f'l <8 J i+ i )i+1 )) 

for s = 1, . . . , mi. In virtue of the lemma we have 

R(R (8) £ m2d )^ +1)i = (8 £m 2 d)^i+ij, (5-9) 

where 

^i+ij+i = diag {exp (-r^ Ji + i,j+i) , exp (-2r^ mi Ji+iJ+i) , • • • , exp (-m 2 r^ mi ^+i,i+l) , 
exp (-r^ J i+ ij +1 ) , exp (-2r^ Ji+ij+i) , . . . , exp (-m 2 r£? mi J i+1 j+i^ 

• ■ • , exp (-rg^ Ji+ij+i) , exp (-2r^ Ji+ij+i) , . . . , exp (-marg^ Ji+ij+i) } . 
Then, taking into account (|5.9p . the relation (|5.8p takes the form 

if + l,i+l^+i,j = (R- 1 E m2d )V} +lij+1 (R E m2d )w} +1>j + 0(r m2 ). (5.10) 
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Next we consider block permutation matrix T = (Tjj), where i = 1, . . . , m 2 and j = 1, . . . , mi, where 
each block Tjj itself consists of blocks Tjj = (Tf. s ), for k = 1, . . . , mi and s = 1, . . . , mi- The blocks 
Tk, s are defined by 

j, _ f #<f, for fe = j, s = i, 
fc ' s \ d , for / j or s / i 

Remark, that the matrix T is in fact a square orthogonal matrix of the order mimid. Also it is worth 

remarking, that 

T (JT 1 ® E m2d ) Vj+u+i (R ® E m2d ) T T = (5.11) 

where 

^}+i,j+l = dia § i ex P (-Hmj ® Ji+ij+i) , exp (-2r7 mi (8) Jj+ij+i) , . . . , exp (-m 2 r7 mi ® Jj+iJ+i)} . 

Then the relation (|5.10p taking into account ()5.1ip becomes 

F} +w w} +ld = 7 Tl b 1 i+1}j+1 7wj +1>j + 0(r m2 ). (5.12) 

Therefore we finished a transformation of Ff +1 j +1 wj +1 j and now let us turn to rK^ +t j +1 wjj +1 in the 
first equation of system (|5.7p . Using (j4.2|) and taking into account (|4.12p we obtain 

rK} +l j+1 w} j+1 = r (E mi <g> 7 m2 ®E d + r^ mi ® £ m2 ® Jj+ij+i) -1 (7^ 8) £ m2 ® Jj+iJ+i) 

= ^m im2 rf + Smi ® 7m 2 ® -tflj+l) ^j+l + O^™ 1 )- 
Again, let i?i be m 2 X m2 matrix, transforming 7 m2 to normal Jordan form 7^ 2 , that is, 
RilmiRi 1 = 7m 2 , where 7 ^ 2 = diag , ^ , . . . , ^ } . 

Then 

ri^ lj+1 ^- +1 = (£ mi ® f^ 1 8) E d ) Ri (E mi ®Ri® E d ) w{ ]+l + 0(h mi ), (5.13) 

where 

#1 = ^E mi m 2 d + ^7mi ® 7m 2 ® J i+l,i+l) ■ 

Remark, that 

= (5.14) 

with .R = diag{^n, R 22 , • • • , flm 2 m 2 }, where 

= ( E mid + (7^ ® JS-ij+i) ) > 

for s = 1, . . . , 7712- Prom (|5,13p and (|5,14p . taking into account lemma U2J we get 

rKl +1J+1 wl j+1 = (E mi ® i^ 1 ® £ d ) T T ^ +lj+1 T(E mi ® ^ ® £ d ) u^ +1 + 0(h m ^), (5.15) 
where 

f ( C 2 ,-1 \ f 2 C 2 ,-1 \ / T -i \ 

Di+ij+i = diag < exp I — J i+1>j+1 I ,exp I — J i+lj+1 I ,...,exp I J mj+1 I , 

f C 2 7 -i \ / ^ 7 -i \ f m ^ 2 7 -i \ 
exp I — J i+1J+1 I ,exp I —J i+W j , ■ ■ ■ , exp I -_J. +1>j+1 I , . . . 

f C \ / 2 C 2 , 1 \ / m lC 2 ,-1 \ \ 

■ • ■ , exp —J i+ i, j+ i h exp — J. +lj . +1 , . . . , exp _ ^ . 
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Remark, that 

(E mi ® R^ 1 ® E d ) 7 T V 2 i+hj+l 7(E mi ® i?i ® E d ) = Df +1)j+1 , (5.16) 

where 

^i+ij+i = dia g jexp (-^7m 2 ® JJ+ij+iJ , exp (-^7m 2 ® ^i,j+i) » • • • > ex P (""^"w ® J m,i+i) } 

Therefore, from (j5. 15[) and (|5.16p we get 

rK} +1)j+l w\ j+l = & i+l!j+l wl j+1 + 0(/i mi ). (5.17) 

To transform the coefficients of the second and third equation of the system (|5.7|) . it is enough to 
make use of (|4.2p and (|4.12p . In a result we obtain 

F? +l>j+1 w? +hj = 0, rKi +ljj+l wl j+1 = w? j+1 + 0{h mi ), 

/f. = + 0(T m >), rK? +hj+1 wl j+1 = 0. (5.18) 

Therefore, a transformation of all coefficients of the system (|5,7p is finished. Let us write down 
transformed difference scheme using (|5.12p . (|5.17p and (|5.18p . We have 

3 

^i+lj'+l = 9i+l,j+l ~ ^ T ^i+lj+l^i+lj ~ ^l+lj+l^ij+l ~ ^ ^+l.? + l(^> T )^\+l,j 

1=1 

- E "l+w T Hi+i + °( hmi ) + °^ m2 ). ( 5 - 19 ) 
/=i 

3 3 

wf+ij+i = ft 2 +i,j+i - Wij+i - E^+i,i+i^' r )^+iJ ~ E ^i+ij+i^' r )^,j'+i + C(/i mi ) ; 

i=i i=i 

3 3 

vfi+ij+i = 9i+i,j+i ~ -J2^t 1J+1 (h,T)w l i+1J -^^ l j+1 (h,T)w l id+1 + 0(T m2 ). 

1=1 1=1 
The system ()5. 19[) represents a canonical form of the difference scheme (|4,2p . In this form one im- 
mediately sees spectral characteristics of matrix coefficients. This information is necessary to prove 
a uniform boundedness of the grid solution in the domain U. Therefore the goal of this section is 
achieved and we are in position to turn to the main part of the paper. 



6. The proof of a stability of the difference scheme 



In this section we prove the following. 

Theorem 6.1. Let in the difference scheme \3. 5|) eigenvalues of the matrices j mi , j m2 and Jj+ij+i 
with given mi, m2 satisfy the condition \4.3\ ) in the domain and j+1 > V i, j and s = 1, . . . , k. 
Let r being a ratio of two steps is a constant. Then the difference scheme \3. 5\) is absolutely stable 
with respect to initial-boundary data and right-hand side. Following estimation for its solution 

\\vi+ij+i\\c(u A ) ^ M i\\fi+i,j+i\\c(u A )+ M 2\\(pi+i\\c(u A ) +M 3 ||V>j+i||c(t/ A ) (6- 1 ) 
is valid, where Mfc are constants and i = 1, . . . , ri\ — 1, j = 1, . . . , n<i — 1. 
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Proof. Remark that the existence of a unique solution of the difference scheme (|3.5p was in fact 
proven in previous section together with its transformation. It remains to prove the second part of 
the theorem. 



The difference scheme (|3,5|) in the previous section was transformed to the form (|5.19|) . To describe 
the structure of difference scheme (|5.19j) . we rewrite it in matrix form 

{L + L(h,r))V = G, (6.2) 

where V is unknown ni^n-dimensional vector, consisting of three block components V = (V 1 V 2 V 3 ), 
where 



V* 



w n u 2i • • • > w m,n 2 



Let us describe in detail system (|6.2p . In (j6.2[l . L is a square block diagonal matrix of the order 
nin^n, that is, L = diagjL 1 , L 2 , L 3 }. In turn each block L k is block two-diagonal matrix L k = (L k j) 
for i = 1, . . . , n\ and j = 1, . . . , r&2- The blocks L k ■ at the main diagonal, also have the block two- 
diagonal form 

/ E s S S . 



L 



'■j 



3^2 E s X 



0, 



o, 

0, 



0., 



(6.3) 



S S S . 

\ Os Os Os . 

Here the parameter s take values d, I and p for k = 1, 2 and 3, respectively. The blocks 3^ of the 
matrix (16.311 have the form 



3~j j — 7 Tl^jl, "Jij — Oi, S'ij — E p . 

The blocks L k situated under the main diagonal, that is, for i = 2, . . . ,rti and j 
following block diagonal form: 



L 



diag 



rv-k 
1 -^1,2) 



, X 



4,ra 2 



} ' X l: 



f) 2 ■ 



X 



Ei. X 



0,. 



(6.4) 
1, have the 

(6.5) 



All the rest blocks L k ■ are zero ones of corresponding sizes. 

The vector G in the system (j6.2|) is as follows: 

G = g + (S + §{h, T ))wo + {Q + Q(h, t))w° + 0(/i mi ) + 0(r ma ), 

with g are nir^n-dimensional vector, consisting of blocks g = (g 1 , g 2 , g 3 ) T of the sizes n\n2d, nin^ 
and niri2P, respectively, 



(6.6) 



Cat ,1) <?1,2; • • • J #1,712' #2,1' #2,2 



~k ~k 



~k 

92,112 



<9ni,l> 9ni,2 



• 9ni 



nil ' 



wq and w° are known vectors of the size n\n,2n the elements of which are defined by initial-boundary 
data f|2.2|) . namely, wq = (wq, w 2 ., Wq) t , where 



and w° = (w ' 1 , w 0,2 , u) 0,3 ) T , where 



w 



1,0' w 2,0i 



■11-2 



w°' k = e n2 ® (w$ tl , w^ 2 ,...,w^ r 
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The matrices S and Q are known square block diagonal ones of the order mr^w, that is, 

5 = diag{S 1 , S 2 , S 3 } and Q = diagjQ 1 , Q 2 , Q 3 } . 

Each square block S k and Q k has the order n\n2d, n\nil and n\n2P, corresponding to the value of k. 
Each block S k is also block diagonal matrix 



diag I'S'i,!, S2 2, ■ ■ ■ , ^n u ni } ' 



where are square blocks of the orders n^d, n^Z and ri2P, respectively k, of the following form: 

S k i = diag {-3*!, S , S ,...,Q S ], 

for i = 1, . . . , m, where s takes values d, I and p, respectively k. Each block Q k is block diagonal 
matrix 

</ = diag<jQ^, S , S ,...,Q,}, 

where 

Qi,i = diag |— 0Ci 1; -3Ci j2 , • • • , — ^i, n2 } • 

Finally, the matrices L(h,r), §(/i,r) and Q(/i,r) in (|6.2p and (|6.6p are quadratic matrices of the 
order ni^n built from the blocks - +1 (h, r) and t). 

To prove the stability property of difference scheme (|3.5p . there is a need to estimate the norm of 
unknown vector V in the system (|6.2p . For this aim we must calculate L . Since L is a block diagonal 
matrix, and each its diagonal block is a block two-diagonal matrix, then we can easily write down the 
matrix L _1 in its explicit form. Each its diagonal block component is of the form 

= A k$k for k = 1,2,3, 



(6.7) 



where A fc are matrices of the order rti^d, nir^Z and n\U2P, respectively value of k, which have the 
form 

/ 



E s 


Os 


Os 


... o s \ 


A k 


E s 




... o s 


A k A k 


AS 




... Os 



A 



vn^=a< n«i=3 a*! ie^a* ... e s ; 

where s take values n2<i, 712/ and n2p, respectively value of k. Every block A k , for % = 2, 
defined as A 4 fc = 1, where 



(6.8) 



, ni is 



/ 



rrk rrk 



0, 

E s 

rrk 

J i,3 



S 

E x 



Os 

Os 

E s j 



(6.9) 



v (-ir 1 ^ 3 ?.. (-i) na+9 iE?=8^ (-ir 2+3 m 2 =4^ si 

where s take values cZ, Z and p, respectively value of k. The matrices <3? fe in (|6.7p have the form 
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Estimate the norms of the matrices (L k )~ 1 . For k = 2 and k = 3, it follows from (|6.4p . (|6.5p and 
l|B77j> - (|B^]l that A* = S and (L fe )" 1 = £ ni „ 2S , where s = d,l for k = 2,3, respectively. Thus, 

||(^)- 1 ||C(C/ A )=1 for k = 2,3. (6.10) 
Let us estimate the norm of the first block component (L 1 ) -1 . It follows from (|6.7p and (|6.8p that 

ll( Ll ) _1 llc(i/ A ) ^ II a1 IIc(c/a)II $1 IIc(c/a)' 

IIAlc^) < 1 + \\tf\\c(U±) + l|Aj||^ A) + • • • + W^W^uly (6-11) 

In virtue of (|6.9|) 

\\^i\\c(U A ) < \\(*lrHc{u A )\\X$Jc<p±)- (6-12) 
In turn for the first multiplier in the right-hand side of the inequality (|6.12p we have 

IK^i) IIc(Ua) - 1+ 11^,^2 Hc(C/ A ) + \\^i,n2^i,n 2 -l\\c{U A ) H H [l^na^r^-l " ' ^,2 llc(C/ A ) " (6-13) 

From (|5.11|) and ()6.4p we find the radius /U of the spectrum of the matrix ■ to get 

M = ma X {|exp (-fc^^) |} (6-14) 
for i = 1, . . . , Tii, j = 1, . . . , n 2 , ki = 1, . . . , TO2, A; 2 = 1, . . . , mi and £3 = 1, . . . , d. We have a pair of 

theorem on spectral decomposition for a power function of matrix, namely, (3~|j) n2 given in (|21j. p. 
155) or theorem from [25], we get 

ll(^rilc(t/A)^°> at n 2 ^oo. (6.15) 
From ()6.15p it follows, that there exists such a value of the power to, for which 

M?hr\\c(u A ) < I, (6-16) 

where m is a suitable constant. Let us consider the product YYj^i ^Ij- With Taylor's formula we 
obtain 

rh+l 

n^ = (^ 2 r + OH- (6-17) 
i=2 

Uf^^lj ■ From dED an d d^TTD it follows that x < 1. Hence, from (joTSj) we 



inequalities: Re(^^ i ) > V k2, which can be verified with the help of Routh-Hurwitz criterion and 
£j 3 . > V k%, which is valid by condition of our theorem. Then, from (|6.14p we get fi < 1. Applying 

rl \ri2 



Denote \ 
get 

where 



C(U A ) 



^)- 1 |b(C/ A ) <ri/(l-x), (6-18) 



v - 1 + ll%llc(c/ A ) + H^IIccc/a) + • • • + II^-IIccc/a)- 



Let us estimate a norm of the matrix %} We learn from |26| that for a constant n x n matrix >1 
spectrum of which lies strictly in the left halfplane it is valid following inequality 

|| exp(L4)|| 2 < v / cexp(-Kt), t > 0, (6.19) 

where k = 1/(2||X||2), c = 1 1 J 5 ^ 1 1 1 2 1 1 1 1 2 and X is Hermitian n x n matrix, being a solution of 
Lyapunov equation XA + A*X = —E n and || • H2 ia a norm of a matrix, agreed with Hermitian norm. 
The constants c and k are clarified in |27|. 
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From the above mentioned assumptions with respect to eigenvalues of the matrix J(x, t) and from 
(|5.16p . (|6.5p and (|6.19|) it follows, that for the matrix %\ ■ it is valid the following estimation: 

||3C^-||c(c/ A ) < Vcexp(-K/r), (6.20) 

where c and k are constants. From (|6.20p it follows that, reducing r, one can achieve a sufficient 
smallness of the elements of Xjj including that the inequality 

-r^-\\X}jc(U.)< 1 - (6-21) 

will be valid. So, from (|6.12p . (|6.18p and (|6.2ip it follows that there are such values of the steps r and 
h for which the inequality ||A|||(7m A ) < 1 holds. Taking into account ([6. lip , we get 

IIA^Io^a) <c, where c = 1/(1 - ||AJ|| c(C7a) ). (6.22) 

Therefore, from (|6.1ip and (|6.18p it follows that 

\\(L 1 r 1 \\c(U A )<c, where c = crj/(l- X ). (6.23) 

These inequalities (|6.10p and (|6.23[) mean uniform boundedness of the matrix in (|6.2p . From (|5.9p 
and (|6.4|) it follows that to estimate a norm of the matrix &j -, we can use the inequality (|6.19p . In a 
result, we get 

\\?}jc(U A ) <^exp(-«r), (6.24) 

where e and R are some constants. Then the matrices S and Q in (|6.6p . in virtue of (|6.2(jp and (|6.24p . 
are also are uniformly boundness at the grid [/a- That is, for sufficiently small values of h and r there 
are such constants p\ and p 2 , for which 

\\S + Hh,T)\\ C (u*)<Pi, \\Q + Q(h,r)\\ c(UA) <p 2 . (6.25) 
From ()6.2p and (I6.6P we obtain 

V = (E nin2h + L- 1 £(/i, r))-^- 1 (g + (S + S(h, r))w + (Q + Q(h, r))w° + 0(/j mi ) + 0(r™ 2 )) . 

(6.26) 

In virtue of a boundness of the matrix L^ 1 , for sufficiently small r and /i, we have 

\\(E ni n 2 h + ^~ 1 ^(^^))" 1 ||c((7a) ^ 

where 77 = 1/(1 — fj) and r] = ||L _1 £(/i, T)\\ C ( UA y Denote 

Mi = TpfjX\\P\\ C{ u A )/(l ~ \\M i+1J+1 \\ C{UA) ), 

M 2 = ppir), M 3 = pp 2 fj, 

where p = max{c, 1} and p±, p 2 as in (|6.25p . Letting h — > and r — )• 0, from (|6.26|) we obtain desired 
estimate (|6.ip . Therefore the theorem is proven. 

Let us remark, the theorem 16.11 is also valid for initial-boundary problems of the form (|2.ip and 
(|2.2p with nondegenerate in a domain of definition matrices A(x,t) and B(x,t). To check this, it is 
enough to put p = and / = 0. 
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7. Numerical experiments 



For numerical solving of boundary problems of the form (|2.ip - (|2.3p with mentioned above condi- 
tions, we have created a program in which the user can set input data, that is, matrix coefficients 
A(x,t), B(x,t), C(x,t) and vector-functions f(x,t), <j)(x) and ip(t); the domain of definition U; the 
values of the steps h and r and the orders of a spline mi and ni2 with respect to each independent 
variable. 

In this section we present the numerical results of solving of some problems of the form ()2. 1 j) - (|2.3[) . 
Remark, that these examples were made only to demonstrate the stability property of the difference 
scheme (13.51). 



Example 7.1. Let us consider the system (|2.ip . in which 



A(x,t) 



V 



1 




1 + xt 







exp(xi) 







\ 







B(x,t) 



( 

exp(sin(i?)) 



xt 



\ 1 



\ 


1 




0/ 



C(x,t) 



f(x,t) 



( exp(tf) \ 
exp(xi)(x exp(xi) + tsin(i?) + 2xt) 

X 

\ 1 ) 

x + t. We know that the exact solution of our system with the data (|7.ip is 

f tf) T . 

xt ^ — 1} matrices 
\ 

















\ 














2xt 




















I? 


1 























1 











V° 














J 



(7.1) 



with $ 

u(x,t) = (x, exp(xi), exp(i?), xt, 1, exp(xi) 
With the help of nondegenerate in compact domain U C {(x,t) G M 2 , 



/ l/exp(xt) 



P(xA) 








V 1 





1 









1/(1 + xt) 







o \ 



1 


°7 



Q 



/o 





1 



1 




-1 












1 



1 

10 

1 



/ 

in which Ji(x,t) 



(7.2) 



the system (|2.ip with the data (|7.ip transforms to the canonical form (f2~^ 
exp(sin(t?))/ exp(xt), J2(x, t) = xt/ (1 + xt) and M = N = O2. To check this, it is enough to multiply 
the system (|2.ip on the left by the matrix P(x, t) and to make a change of variable: v(x, t) = Qu(x, t). 
Remark, that in arbitrary compact domain U the first condition of the theorem 12.21 is not fulfilled. To 
check this, it is enough to write down the characteristic equation of the system (|2.ip with the data 
(|7.ip and to see that its roots have common values in U. Nevertheless, the system (|2,ip with the data 
(|7.1|) can be cast to the form (|2.4|) and therefore the conditions of the theorem 16. II are fulfilled. Thus, 
in this case the difference scheme (|3.5p is stable in the domain U. We show the results of numerical 
solving in the table 1. As stability estimate we take in this case the value of absolute error of the 
solution Au = \\u itj - v itj \\ C ( UA y 
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Table 1. 



N 


h 


r 


to 


T 


x 


X 


mi 


777-2 


Au 


1 


io- 1 


io- 1 





1 





1 


2 


2 


2.07 x 10~ 2 


2 


io- 1 


io- 1 





1 





1 


3 


2 


2.07 x 10~ 2 


4 


io- 1 


io- 1 





1 





1 


2 


3 


1.96 x 10~ 3 


5 


io- 1 


io- 1 





1 





1 


3 


3 


1.96 x 10~ 3 


6 


io- 1 


io- 1 





1 





1 


4 


3 


1.96 x 10~ 3 


7 


io- 1 


io- 1 





1 





1 


3 


4 


1.91 x IO" 4 


8 


io- 1 


io- 1 





1 





1 


4 


4 


1.91 x IO" 4 


9 


io- 1 


io- 1 





1 





1 


5 


5 


1.91 x 10~ 5 


10 


io- 1 


io- 1 





1 





1 


6 


6 


1.91 x 10~ 6 


11 


io- 1 


io- 1 





1 





1 


7 


7 


1.93 x 10~ 7 


12 


io- 2 


io- 1 





1 





1 


2 


2 


2.07 x 10~ 2 


13 


io- 1 


io- 2 





1 





1 


2 


2 


1.28 x 10~ 3 


14 


io- 2 


io- 2 





1 





1 


2 


2 


1.96 x IO" 4 


15 


5 x 10~ 3 


5 x IO" 3 





1 





1 


2 


2 


4.95 x 10~ 5 


16 


io- 1 


KT 1 





1 





2 


2 


2 


7.14 x 10~ 2 


17 


io- 1 


io- 1 





1 





2 


3 


3 


1.31 x 10~ 2 


18 


io- 1 


io- 1 





1 





2 


4 


4 


2.52 x 10~ 3 


19 


IO' 1 


io- 1 





1 





2 


5 


5 


4.99 x IO" 4 



From the table 1 we see, that to achieve more accuracy of numerical solving, in some cases, it is 
enough to increase a degree of the spline with respect to one variable. This is important, because a 
degree of the spline impact on the dimension of the difference scheme (|3.5p and therefore on a speed 
of calculations. What we learn from the test 17 is that with using third-degree spline the accuracy of 
the numerical solutions in the domain U = [0; 1] x [0; 2] is of the second degree. This is explained by 
the impact of a large Lipchitz constant. 

Example 7.2. Let us consider the system (|2.ip with multiple characteristic curves with the following 
data 

A(x,t) = diag{£ 5 , 0, 1}, B(x,t)=diag{Ji(x,t), J 2 {x,t), 1, 0}, 



exp(tf) 1 
Ji(x,t) = | exp(tf) 1 

exp(tf) 



1 + t exp(x) 1 

l + iexp(z) 



C(x,t) 



( 


x 2 + t 





1 


l + xt 


— exp(i9) 





\ 







;r 2 


xt 








1 


7? 




1 











1 


xt 













1 







































x exp($) 


1 


xH 





i) 








\ 








exp(i?) 








1 


/ 



(7.3) 
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/ exp(2tf) + (x 2 + t) exp(#) + 2xt + (1 + xt)(x - t) + 1 \ 
exp(i?) + 2t + x 3 + x 2 t + 3x 2 i 2 + x exp(t) + x 3 t + 1 
(2i + 1) exp(tf) + 2x + x 2 t exp(t) + 1 
f(x,t)= 2xt + texp(x) 


x exp(2??) + exp(t) + 2?? + 2x 3 t 2 
^ x 2 + 2xi exp(i9) + x exp(i) y 

Here $ is defined in example 1. The exact solution of (|2.1|) with the data (|7.3|) is 

u(x, t) = (exp(-#), i?, 2xt, x — t, 1, xexp(i), x 2 t) T . 

This system is given in canonical form. It has multiple nontrivial characteristic curves 

Ai,2,3 = -exp(-#), A 4j5 = -1/(1 + texp(x)) 

and one zero characteristics Xq = 0. In an arbitrary compact domain U C U = {(x,t), texp(x) > 0}, 
the conditions of theorems 12. 2l and l6.1l and remark [2T3l are valid. Therefore, the difference scheme (|3.5|) 
is stable in any domain U. The results of numerical solving are given in the table 2. 

Table 2. 



N 


h 


T 


to 


T 


x 


X 


mi 


m 2 


An 


1 


io- 1 


10" 


1 





1 





1 


2 


2 


3.54 x 1Q-' Z 


2 


10" 1 


10" 


1 





1 





1 


3 


2 


4.96 x 10~ 3 


3 


io- 1 


io- 


1 





1 





1 


2 


3 


3.04 x 10~ 2 


4 


io- 1 


10" 


1 





1 





1 


3 


3 


3.34 x 10~ 3 


5 


10" 1 


io- 


1 





1 





1 


4 


3 


4.78 x 10~ 4 


6 


io- 1 


io- 


1 





1 





1 


3 


4 


2.91 x 10~ 3 


7 


10" 1 


io- 


1 





1 





1 


4 


4 


3.23 x 10~ 4 


8 


io- 1 


io- 


1 





1 





1 


5 


3 


6.19 x 10~ 4 


9 


10" 1 


io- 


1 





1 





1 


5 


4 


4.76 x IO" 5 


10 


10" 1 


io- 


1 





1 





1 


5 


5 


3.21 x 10~ 5 


11 


10" 2 


io- 


2 





1 





1 


2 


2 


3.23 x 10~ 4 


12 


io- 1 


io- 


1 





2 





2 


2 


2 


7.97 


13 


io- 1 


io- 


1 





2 





2 


3 


3 


7.64 x IO" 1 


14 


io- 1 


10" 


1 





2 





2 


4 


4 


7.46 x IO" 2 



Remark, that in the tests 12, 13 and 14 on the error has influenced a large Lipchitz constant. 



8. Conclusion 

Numerous actual calculations show that spline collocation difference scheme given in this paper is 
a quite effective and gives sufficient accuracy. But, as was mentioned above, the difference scheme 
(|3.5p in fact being a linear system of algebraic equations has the order m\m2n. Therefore increase 
of the degree of a used spline yields an increase of an order of this system. Since our system by 
assumption has general form then to solve it one needs to use, generally speaking, universal methods 
and therefore procedure of numerical solving might be quite laborious. In this case there is a need 
of a parallelization of numerical calculations. Following remark is in order. The point is that many 
actual problems of mathematical physics and mechanics are described by systems of partial differential 
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algebraic equations, involving more than two independent variables. Examples are given by Sobolev's 
system, linearized Navier-Stokes system etc. [Ij. Therefore the method given in this work could be 
more marketable if one generalizes it to the case of many independent variables and this is the problem 
to be addressed. 
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